Nigerian-TCGA Differential Expression Analysis - Exploration of solutions to batch effect normalization
Padma Sheila Rajagopal, MD MPH
7/24/2019
Last updated: 2020-06-07
Checks: 5 2
Knit directory: ~/Research-Local/2019-rnaseq/TCGA-Nigerian-RNAseq/
This reproducible R Markdown analysis was created with workflowr (version 1.4.0). The Checks tab describes the reproducibility checks that were applied when the results were created. The Past versions tab lists the development history.
The R Markdown file has unstaged changes. To know which version of the R Markdown file created these results, you’ll want to first commit it to the Git repo. If you’re still working on the analysis, you can ignore this warning. When you’re finished, you can run wflow_publish to commit the R Markdown file and build the HTML.
Great job! The global environment was empty. Objects defined in the global environment can affect the analysis in your R Markdown file in unknown ways. For reproduciblity it’s best to always run the code in an empty environment.
The command set.seed(20190113) was run prior to running the code in the R Markdown file. Setting a seed ensures that any results that rely on randomness, e.g. subsampling or permutations, are reproducible.
Great job! Recording the operating system, R version, and package versions is critical for reproducibility.
- DE test: Nigerian basal vs Nigerian Her2
- DE test: Nigerian basal vs TCGA white basal and TCGA black basal for RUV
- DE test: Nigerian basal vs TCGA white basal
- DE test: TCGA black and TCGA white for RUV
- DE test: TCGA black and TCGA white
- Differential expression analysis -> batch effect comparison
- Hand code value check of Her2
- Methods cross check
- Nigerian group cross check
- RUVSeq use of housekeeping genes for normalization
To ensure reproducibility of the results, delete the cache directory NigerianTCGArawcountsDeSeq2-pc2_cache and re-run the analysis. To have workflowr automatically delete the cache directory prior to building the file, set delete_cache = TRUE when running wflow_build() or wflow_publish().
Great job! Using relative paths to the files within your workflowr project makes it easier to run your code on other machines.
Great! You are using Git for version control. Tracking code development and connecting the code version to the results is critical for reproducibility. The version displayed above was the version of the Git repository at the time these results were generated.
Note that you need to be careful to ensure that all relevant files for the analysis have been committed to Git prior to generating the results (you can use wflow_publish or wflow_git_commit). workflowr only checks the R Markdown file, but you know if there are other scripts or data files that it depends on. Below is the status of the Git repository when the results were generated:
Ignored files:
Ignored: .DS_Store
Ignored: .Rproj.user/
Ignored: analysis/.DS_Store
Ignored: analysis/.Rhistory
Ignored: code/.DS_Store
Ignored: docs/.DS_Store
Ignored: docs/figure/.DS_Store
Ignored: docs/figure/NigerianTCGArawcountsDeSeq2-pc2.Rmd/.DS_Store
Ignored: plots/.DS_Store
Untracked files:
Untracked: NigerianTCGArawcountsDeSeq2-proteincoding-IHC.Rmd
Untracked: NigerianTCGArawcountslimma-voomDE-PAM50-linc.Rmd
Untracked: NigerianTCGArawcountslimma-voomDE-PAM50.Rmd
Untracked: _site.yml
Untracked: about.Rmd
Untracked: analysis/NigerianTCGArawcountsDeSeq2-pc2_cache/
Untracked: analysis/figure/
Untracked: docs/figure/NigerianTCGArawcountslimma-voomDE-PAM50.Rmd/
Untracked: index.Rmd
Untracked: plots/Nanostring-genes-in-RNAseq-results.jpeg
Unstaged changes:
Modified: .Rhistory
Modified: README.md
Modified: analysis/NigerianTCGArawcountsDeSeq2-pc2.Rmd
Deleted: analysis/NigerianTCGArawcountsDeSeq2-proteincoding-IHC.Rmd
Deleted: analysis/NigerianTCGArawcountslimma-voomDE-PAM50.Rmd
Deleted: analysis/_site.yml
Deleted: analysis/about.Rmd
Deleted: analysis/index.Rmd
Modified: code/proteincodingparse.Rmd
Note that any generated files, e.g. HTML, png, CSS, etc., are not included in this status report because it is ok for generated content to have uncommitted changes.
These are the previous versions of the R Markdown and HTML files. If you’ve configured a remote Git repository (see ?wflow_git_remote), click on the hyperlinks in the table below to view them.
| File | Version | Author | Date | Message |
|---|---|---|---|---|
| Rmd | 83061c6 | Sheila Rajagopal | 2019-08-19 | Differential expression with IHC and LOH subanalysis |
| Rmd | 0a03483 | Sheila Rajagopal | 2019-06-19 | Gene check for all normalization methods |
| Rmd | c56450c | Sheila Rajagopal | 2019-06-18 | updated RUVseq PCA |
| Rmd | c7160b4 | Sheila Rajagopal | 2019-06-17 | Updated quantile normalization |
| Rmd | 1b2b589 | Sheila Rajagopal | 2019-05-08 | batch effect assessment |
#Translation from HTSeq raw counts -> DESeq2 table of raw counts I have 84 TCGA patients (70 white or black) with whole-genome sequencing data and RNAseq data as well as 96 Nigerian patients with RNA-seq data. Raw counts were initially processed using HTSeq, so HTSeq raw counts (with one file per patient) are being formatted for use with DESeq2.
sampleConditionPAM50
sampleConditionrace Basal Her2 LumA LumB Normal PAM_other
Nigerian 41 27 14 11 3 0
TCGA_black 23 0 4 4 0 0
TCGA_other 0 0 0 0 0 14
TCGA_white 17 5 8 9 0 0
Warning: The above code chunk cached its results, but it won’t be re-run if previous chunks it depends on are updated. If you need to use caching, it is highly recommended to also set knitr::opts_chunk$set(autodep = TRUE) at the top of the file (in a chunk that is not cached). Alternatively, you can customize the option dependson for each individual chunk that is cached. Using either autodep or dependson will remove this warning. See the knitr cache options for more details.
Warning: The above code chunk cached its results, but it won’t be re-run if previous chunks it depends on are updated. If you need to use caching, it is highly recommended to also set knitr::opts_chunk$set(autodep = TRUE) at the top of the file (in a chunk that is not cached). Alternatively, you can customize the option dependson for each individual chunk that is cached. Using either autodep or dependson will remove this warning. See the knitr cache options for more details.
Figure 1: Graphical representation of the variance-stabilizing transformation relative to the standard transformation of log2(x+1). This reduces variance purely due to low read counts. VST is effectively employed in DESeq2. This transformation will be the input for subsequent visualizations (ONLY VISUALIZATIONS) in this section.
Warning: The above code chunk cached its results, but it won’t be re-run if previous chunks it depends on are updated. If you need to use caching, it is highly recommended to also set knitr::opts_chunk$set(autodep = TRUE) at the top of the file (in a chunk that is not cached). Alternatively, you can customize the option dependson for each individual chunk that is cached. Using either autodep or dependson will remove this warning. See the knitr cache options for more details.
Figure 2: Plot of mean to variance to confirm that variance shrinking does work. This is another visualization of the approach noted in Figure 1.
Warning: The above code chunk cached its results, but it won’t be re-run if previous chunks it depends on are updated. If you need to use caching, it is highly recommended to also set knitr::opts_chunk$set(autodep = TRUE) at the top of the file (in a chunk that is not cached). Alternatively, you can customize the option dependson for each individual chunk that is cached. Using either autodep or dependson will remove this warning. See the knitr cache options for more details.
PCA 1-1:PCA demonstrating the raw difference between Nigerian and TCGA patients with within-batch demonstration. Nigerian patients do not appear to have significant differneces based on their batches. (TCGA patients are known not to have significant difference based on batch performed.)
There is a factor separating the TCGA patients into 2 groups that is not wholly ancestry dependent.
Of note, the PCA as performed by only has 33% variance explained with two variables when transformation is performed.
Warning: The above code chunk cached its results, but it won’t be re-run if previous chunks it depends on are updated. If you need to use caching, it is highly recommended to also set knitr::opts_chunk$set(autodep = TRUE) at the top of the file (in a chunk that is not cached). Alternatively, you can customize the option dependson for each individual chunk that is cached. Using either autodep or dependson will remove this warning. See the knitr cache options for more details.
PCA 1-2:PCA demonstrating association with breast cancer subtype and gene expression. This is noted in both groups and consistent with known biology. Again, only 32% variance is explained with 2 principal components.
Warning: The above code chunk cached its results, but it won’t be re-run if previous chunks it depends on are updated. If you need to use caching, it is highly recommended to also set knitr::opts_chunk$set(autodep = TRUE) at the top of the file (in a chunk that is not cached). Alternatively, you can customize the option dependson for each individual chunk that is cached. Using either autodep or dependson will remove this warning. See the knitr cache options for more details.
PCA 1-3: Integrated PCA demonstrating combined effect of breast cancer subtype and race. Subtype appears to account for the clustering seen in the TCGA group. Again, only 32% variance is explained with 2 principal components (with top 1000 genes).
Warning: The above code chunk cached its results, but it won’t be re-run if previous chunks it depends on are updated. If you need to use caching, it is highly recommended to also set knitr::opts_chunk$set(autodep = TRUE) at the top of the file (in a chunk that is not cached). Alternatively, you can customize the option dependson for each individual chunk that is cached. Using either autodep or dependson will remove this warning. See the knitr cache options for more details.
Heatmap 1-1: This heatmap is only using the top 250 genes and is competely heirarchically clustered. Results from this are consistent with the PCA.
#Specific issues associated with our batches -We do not know the extent to which general RNA expression varies within the Nigerian population compared to within TCGA populations because of ethnicity / common ancestral variants -Biological differences are potentially confounded by batch effect (as Nigerians are a separate batch from TCGA); although no evidence of within-batch effect
We will discuss the following methods to manage the batch effect and our reasoning for our chosen approach. (1) sva (1 principal component) (2) RUVseq of housekeeping genes (k=1, with 2 and 6 housekeeping genes) (3) Quantile normalization (4) Limma removeBatchEffect
#Method 1: sva to account for batch effect related to the Nigerian/TCGA samples We are using SVA to create an unsupervised hypothetical “batch effect” variable that is then used to account for testing differences between the two groups.
Number of significant surrogate variables is: 1
Iteration (out of 5 ):1 2 3 4 5
Warning: The above code chunk cached its results, but it won’t be re-run if previous chunks it depends on are updated. If you need to use caching, it is highly recommended to also set knitr::opts_chunk$set(autodep = TRUE) at the top of the file (in a chunk that is not cached). Alternatively, you can customize the option dependson for each individual chunk that is cached. Using either autodep or dependson will remove this warning. See the knitr cache options for more details.
Use of sva to create an unsupervised “batch effect” variable that effectively has full interaction with the Nigerian vs. TCGA variable does cause issues – performing differential expression analysis under these conditions does lead to not converging in beta.
I am performing an intial differential expression analysis with a batch-uncorrected version of the counts as well as the version that incorporates SVA-based batch correction.
[1] "Intercept"
[2] "sampleCondition_Nigerian.Her2_vs_Nigerian.Basal"
[3] "sampleCondition_Nigerian.LumA_vs_Nigerian.Basal"
[4] "sampleCondition_Nigerian.LumB_vs_Nigerian.Basal"
[5] "sampleCondition_Nigerian.Normal_vs_Nigerian.Basal"
[6] "sampleCondition_TCGA_black.Basal_vs_Nigerian.Basal"
[7] "sampleCondition_TCGA_black.LumA_vs_Nigerian.Basal"
[8] "sampleCondition_TCGA_black.LumB_vs_Nigerian.Basal"
[9] "sampleCondition_TCGA_other.PAM_other_vs_Nigerian.Basal"
[10] "sampleCondition_TCGA_white.Basal_vs_Nigerian.Basal"
[11] "sampleCondition_TCGA_white.Her2_vs_Nigerian.Basal"
[12] "sampleCondition_TCGA_white.LumA_vs_Nigerian.Basal"
[13] "sampleCondition_TCGA_white.LumB_vs_Nigerian.Basal" [1] "Intercept"
[2] "be1"
[3] "sampleCondition_Nigerian.Her2_vs_Nigerian.Basal"
[4] "sampleCondition_Nigerian.LumA_vs_Nigerian.Basal"
[5] "sampleCondition_Nigerian.LumB_vs_Nigerian.Basal"
[6] "sampleCondition_Nigerian.Normal_vs_Nigerian.Basal"
[7] "sampleCondition_TCGA_black.Basal_vs_Nigerian.Basal"
[8] "sampleCondition_TCGA_black.LumA_vs_Nigerian.Basal"
[9] "sampleCondition_TCGA_black.LumB_vs_Nigerian.Basal"
[10] "sampleCondition_TCGA_other.PAM_other_vs_Nigerian.Basal"
[11] "sampleCondition_TCGA_white.Basal_vs_Nigerian.Basal"
[12] "sampleCondition_TCGA_white.Her2_vs_Nigerian.Basal"
[13] "sampleCondition_TCGA_white.LumA_vs_Nigerian.Basal"
[14] "sampleCondition_TCGA_white.LumB_vs_Nigerian.Basal"
Warning: The above code chunk cached its results, but it won’t be re-run if previous chunks it depends on are updated. If you need to use caching, it is highly recommended to also set knitr::opts_chunk$set(autodep = TRUE) at the top of the file (in a chunk that is not cached). Alternatively, you can customize the option dependson for each individual chunk that is cached. Using either autodep or dependson will remove this warning. See the knitr cache options for more details.
##sva -> differential expression analysis overall results
MA Plot 1: Original differential expression
MA Plot 2: Batch-corrected differential expression
Warning: The above code chunk cached its results, but it won’t be re-run if previous chunks it depends on are updated. If you need to use caching, it is highly recommended to also set knitr::opts_chunk$set(autodep = TRUE) at the top of the file (in a chunk that is not cached). Alternatively, you can customize the option dependson for each individual chunk that is cached. Using either autodep or dependson will remove this warning. See the knitr cache options for more details.
An MAPlot demonstrates log2 fold changes attributable to the condition over the mean of normalized counts for all the samples in the DESeqDataSet. Points will be colored red if the adjusted p value is less than 0.1. Points which fall out of the window are plotted as open triangles pointing either up or down.
We see an overrepresentation of significant negative log fold changes of expression in the Nigerian patients relative to the TCGA white patients, but this is exacerbated after the sva correction.
##sva -> Inter-TCGA test
[1] 19724[1] 336
[1] 19724[1] 536
Warning: The above code chunk cached its results, but it won’t be re-run if previous chunks it depends on are updated. If you need to use caching, it is highly recommended to also set knitr::opts_chunk$set(autodep = TRUE) at the top of the file (in a chunk that is not cached). Alternatively, you can customize the option dependson for each individual chunk that is cached. Using either autodep or dependson will remove this warning. See the knitr cache options for more details.
The histogram of unadjusted p-values for both uncorrected and SVA corrected appears anticonservative. Both volcano plots demonstrate a tendency for increased log fold change in TCGA white relative to TCGA black, but the batch corrected volcano plot shows greater dispersion / more significance in findings that is concerning for inappropriate exaggeration of effect of positive differential expression.
The following is a demonstration example of comparing 32 Nigerian basal cases and 17 TCGA white basal cases.
##sva -> Basal Nigerian-TCGA white (cross population) test
[1] 19724[1] 8196
[1] 19724[1] 8372
Warning: The above code chunk cached its results, but it won’t be re-run if previous chunks it depends on are updated. If you need to use caching, it is highly recommended to also set knitr::opts_chunk$set(autodep = TRUE) at the top of the file (in a chunk that is not cached). Alternatively, you can customize the option dependson for each individual chunk that is cached. Using either autodep or dependson will remove this warning. See the knitr cache options for more details.
Although unadjusted p-value histogram appears anti-conservative again for both non batch corrected and sva corrected, volcano plot with batch correction appears to have an exaggered effect with numerous genes demonstrating >50x log fold change with significance. The volcano plot without batch correction is more consistent with expected findings, making the sva method concerning for use.
#Method 2: RUVseq of housekeeping genes Based on the paper by Krasnov et al. (https://www.frontiersin.org/articles/10.3389/fgene.2019.00097/full), we examined the housekeeping genes UBC and PUM1 across the Nigerian and TCGA samples to assess their differential expression. We also used the paper by Tilli et al. (https://bmcgenomics.biomedcentral.com/articles/10.1186/s12864-016-2946-1) to assess the efficacy of adding CCSER2, SYMPK and ANKRD17 for housekeeping comparison.
UBC: ENSG00000150991.10 PUM1: ENSG00000134644.11 CCSER2: ENSG00000107771.11 SYMPK: ENSG00000125755.14 ANKRD17: ENSG00000132466.13 -> Given the broad distribution of expression of this gene, this did not appear to be a good standard for comparison in our cohort.
log2 fold change (MLE): sampleCondition Nigerian.Basal vs TCGA_white.Basal
Wald test p-value: sampleCondition Nigerian.Basal vs TCGA_white.Basal
DataFrame with 1 row and 6 columns
baseMean log2FoldChange lfcSE
<numeric> <numeric> <numeric>
ENSG00000150991.10 6998.82317304847 1.16987416588514 0.129891270044962
stat pvalue padj
<numeric> <numeric> <numeric>
ENSG00000150991.10 9.00656499455415 2.12611448406579e-19 6.70281667311056e-18log2 fold change (MLE): sampleCondition Nigerian.Basal vs TCGA_white.Basal
Wald test p-value: sampleCondition Nigerian.Basal vs TCGA_white.Basal
DataFrame with 1 row and 6 columns
baseMean log2FoldChange lfcSE
<numeric> <numeric> <numeric>
ENSG00000134644.11 1774.30511432482 -0.671881253019561 0.109756080788103
stat pvalue padj
<numeric> <numeric> <numeric>
ENSG00000134644.11 -6.12158568523149 9.26486734899955e-10 7.74755265931127e-09log2 fold change (MLE): sampleCondition Nigerian.Basal vs TCGA_white.Basal
Wald test p-value: sampleCondition Nigerian.Basal vs TCGA_white.Basal
DataFrame with 1 row and 6 columns
baseMean log2FoldChange lfcSE
<numeric> <numeric> <numeric>
ENSG00000107771.11 666.733763208679 -0.611348030981247 0.134741414189675
stat pvalue padj
<numeric> <numeric> <numeric>
ENSG00000107771.11 -4.53719470481923 5.70074702127162e-06 2.4568906667344e-05log2 fold change (MLE): sampleCondition Nigerian.Basal vs TCGA_white.Basal
Wald test p-value: sampleCondition Nigerian.Basal vs TCGA_white.Basal
DataFrame with 1 row and 6 columns
baseMean log2FoldChange lfcSE
<numeric> <numeric> <numeric>
ENSG00000125755.14 1792.51275479681 0.185943941133008 0.17491360446888
stat pvalue padj
<numeric> <numeric> <numeric>
ENSG00000125755.14 1.06306162803986 0.287754004379429 0.381899263898067
Warning: The above code chunk cached its results, but it won’t be re-run if previous chunks it depends on are updated. If you need to use caching, it is highly recommended to also set knitr::opts_chunk$set(autodep = TRUE) at the top of the file (in a chunk that is not cached). Alternatively, you can customize the option dependson for each individual chunk that is cached. Using either autodep or dependson will remove this warning. See the knitr cache options for more details.
UBC and PUM1 (as well as CCSER2) demonstrate signifcantly different distribution between the Nigerian and TCGA samples (although the log2 fold change is less for the PUM1 gene). SYMPK does not have significantly different distribution. This suggests that normalization of the housekeeping genes can be used as a means of comparing the RNA read counts.
Warning: The above code chunk cached its results, but it won’t be re-run if previous chunks it depends on are updated. If you need to use caching, it is highly recommended to also set knitr::opts_chunk$set(autodep = TRUE) at the top of the file (in a chunk that is not cached). Alternatively, you can customize the option dependson for each individual chunk that is cached. Using either autodep or dependson will remove this warning. See the knitr cache options for more details.
The RUV normalization demonstrates a different PCA with some clustering based on Nigerian vs. TCGA and some loss of clustering based on subtype relative to VSD normalization or SVA batch effect correction.
##RUVseq -> differential expression analysis overall results
MA Plot: RUVSeq-corrected differential expression
Warning: The above code chunk cached its results, but it won’t be re-run if previous chunks it depends on are updated. If you need to use caching, it is highly recommended to also set knitr::opts_chunk$set(autodep = TRUE) at the top of the file (in a chunk that is not cached). Alternatively, you can customize the option dependson for each individual chunk that is cached. Using either autodep or dependson will remove this warning. See the knitr cache options for more details.
Again, we see an overrepresentation of significant negative log fold changes of expression in the Nigerian patients relative to the TCGA white patients that is exacerbated by the RUV correction.
##RUVseq -> Inter-TCGA test
[1] 19724[1] 500
Warning: The above code chunk cached its results, but it won’t be re-run if previous chunks it depends on are updated. If you need to use caching, it is highly recommended to also set knitr::opts_chunk$set(autodep = TRUE) at the top of the file (in a chunk that is not cached). Alternatively, you can customize the option dependson for each individual chunk that is cached. Using either autodep or dependson will remove this warning. See the knitr cache options for more details.
Again, we see an overrepresentation of differentially expressed genes on the positive side relative to the negative, concerning for the same issue as seen in SVA batch correction.
##RUVSeq -> Basal Nigerian-TCGA white and TCGA black (cross population) tests
[1] 19724[1] 7339
[1] 19724[1] 6506
Warning: The above code chunk cached its results, but it won’t be re-run if previous chunks it depends on are updated. If you need to use caching, it is highly recommended to also set knitr::opts_chunk$set(autodep = TRUE) at the top of the file (in a chunk that is not cached). Alternatively, you can customize the option dependson for each individual chunk that is cached. Using either autodep or dependson will remove this warning. See the knitr cache options for more details.
Although unadjusted p-value histogram appears anti-conservative again for RUV data, volcano plot with RUV correction appears to have the same exaggered effect as seen with the SVA batch effect.
There is some reassurance that both methods are finding genes in common.
#Method 3: Quantile normalization Quantile normalization is another technique that can be used to make distributions comparable. To quantile normalize two or more distributions to each other, without a reference distribution, one must sort the values within the distribution, then set to the average (usually, arithmetic mean) of the distributions. So the highest value in all cases becomes the mean of the highest values, the second highest value becomes the mean of the second highest values, and so on.
Warning: The above code chunk cached its results, but it won’t be re-run if previous chunks it depends on are updated. If you need to use caching, it is highly recommended to also set knitr::opts_chunk$set(autodep = TRUE) at the top of the file (in a chunk that is not cached). Alternatively, you can customize the option dependson for each individual chunk that is cached. Using either autodep or dependson will remove this warning. See the knitr cache options for more details.
The quantile normalization demonstrates a PCA that has similar clustering and % explanations relative to VSD normalization, but the differences between groups have been narrowed.
##Quantile normalization -> DE analysis and Inter-TCGA test
TCGAblackvsTCGAwhite
Down 0
NotSig 18969
Up 1
Warning: The above code chunk cached its results, but it won’t be re-run if previous chunks it depends on are updated. If you need to use caching, it is highly recommended to also set knitr::opts_chunk$set(autodep = TRUE) at the top of the file (in a chunk that is not cached). Alternatively, you can customize the option dependson for each individual chunk that is cached. Using either autodep or dependson will remove this warning. See the knitr cache options for more details.
It has been previously demonstrated that there are expression difference between ethnic groups within the TCGA cohort, and this technique appears to drastically minimize this effect.
##Quantile normalization -> Basal Nigerian-TCGA white (cross population) test
NigerianvsTCGA
Down 3734
NotSig 12329
Up 3364
Warning: The above code chunk cached its results, but it won’t be re-run if previous chunks it depends on are updated. If you need to use caching, it is highly recommended to also set knitr::opts_chunk$set(autodep = TRUE) at the top of the file (in a chunk that is not cached). Alternatively, you can customize the option dependson for each individual chunk that is cached. Using either autodep or dependson will remove this warning. See the knitr cache options for more details.
As per Wei Shi’s preferred method: Counts were converted to log2 counts per million, quantile normalized and precision weighted with the ‘voom’ function of the limma package46,47. A linear model was fitted to each gene, and empirical Bayes moderated t-statistics were used to assess differences in expression48. Empirical Bayes moderated-t P values were computed relative to a fold-change cutoff of 1.5-fold. P values were adjusted to control the global FDR across all comparisons with the ‘global’ option of the limma package. Genes were considered differentially expressed if they had an FDR of 0.05.
##Limma’s removeBatchEffect for visualization only
Warning: The above code chunk cached its results, but it won’t be re-run if previous chunks it depends on are updated. If you need to use caching, it is highly recommended to also set knitr::opts_chunk$set(autodep = TRUE) at the top of the file (in a chunk that is not cached). Alternatively, you can customize the option dependson for each individual chunk that is cached. Using either autodep or dependson will remove this warning. See the knitr cache options for more details.
Warning: The above code chunk cached its results, but it won’t be re-run if previous chunks it depends on are updated. If you need to use caching, it is highly recommended to also set knitr::opts_chunk$set(autodep = TRUE) at the top of the file (in a chunk that is not cached). Alternatively, you can customize the option dependson for each individual chunk that is cached. Using either autodep or dependson will remove this warning. See the knitr cache options for more details.
This heatmap after batch correction with limma’s removeBatchEffect method is using the top 250 genes and is competely heirarchically clustered. These clusters are difficult to separate from each other, and it is not clear that the Limma removeBatchEffect method does not overcompensate and competely flatten existing biological signal.
In addition, the voom method is considered more appropriate for data with significant noise / abnormalities.
#Internal data validity assessments We will test the internal validity of using DESeq2 using comparisons between subtypes, comparisons within the Nigerian cohort. This will not use the batch effect corrected approach with SVA given the demonstrations we showed earlier of the inappropriate exaggeration of positive findings.
##Inter-subtype gene expression testing
[1] 19724[1] 3545
Warning: The above code chunk cached its results, but it won’t be re-run if previous chunks it depends on are updated. If you need to use caching, it is highly recommended to also set knitr::opts_chunk$set(autodep = TRUE) at the top of the file (in a chunk that is not cached). Alternatively, you can customize the option dependson for each individual chunk that is cached. Using either autodep or dependson will remove this warning. See the knitr cache options for more details.
The unadjusted p-value histogram and volcano plot appear appropriate without batch correction in the cross-subtype test. It is reassuring that we are seeing more normal volcano plots on the within-assay comparison.
log2 fold change (MLE): sampleCondition Nigerian.Basal vs Nigerian.Her2
Wald test p-value: sampleCondition Nigerian.Basal vs Nigerian.Her2
DataFrame with 1 row and 6 columns
baseMean log2FoldChange lfcSE
<numeric> <numeric> <numeric>
ENSG00000141736.9 26255.4348341208 -3.31609248417085 0.337356712337551
stat pvalue padj
<numeric> <numeric> <numeric>
ENSG00000141736.9 -9.82963244215175 8.39245146092649e-23 7.43361388151564e-20
Welch Two Sample t-test
data: tempassaybasal and tempassayher2
t = -4.0023, df = 34.362, p-value = 0.0003174
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
-287127.31 -93786.05
sample estimates:
mean of x mean of y
27782.88 218239.56
Warning: The above code chunk cached its results, but it won’t be re-run if previous chunks it depends on are updated. If you need to use caching, it is highly recommended to also set knitr::opts_chunk$set(autodep = TRUE) at the top of the file (in a chunk that is not cached). Alternatively, you can customize the option dependson for each individual chunk that is cached. Using either autodep or dependson will remove this warning. See the knitr cache options for more details.
The cross subtype comparison did find a statistically significant difference between HER2 expression between the Nigerian basal and Her tumors, which is appropriate when the data is visualized as well.
We will now assess arbitrary within-Nigerian samples to identify any other within sample sources of bias.
##Arbitrary grouping of Nigerian patients cross-checked with each other
sampleConditionPAM50
batchval Basal Her2 LumA LumB Normal PAM_other
batch1 10 6 1 2 1 0
batch23 11 12 6 5 1 0
batch4 7 3 4 0 0 0
batch5 13 6 3 4 1 0
batchT 40 5 12 13 0 14[1] 2265
Warning: The above code chunk cached its results, but it won’t be re-run if previous chunks it depends on are updated. If you need to use caching, it is highly recommended to also set knitr::opts_chunk$set(autodep = TRUE) at the top of the file (in a chunk that is not cached). Alternatively, you can customize the option dependson for each individual chunk that is cached. Using either autodep or dependson will remove this warning. See the knitr cache options for more details.
#Methods cross check To compare results between methods to identify the best approach to normalization, we noted the top and bottom log2foldchange genes per method.
[1] "SVA Batch effect correction"[1] "RUV housekeeping gene Batch effect correction"[1] "Limma-voom with quantile normalization batch effect correction"
Warning: The above code chunk cached its results, but it won’t be re-run if previous chunks it depends on are updated. If you need to use caching, it is highly recommended to also set knitr::opts_chunk$set(autodep = TRUE) at the top of the file (in a chunk that is not cached). Alternatively, you can customize the option dependson for each individual chunk that is cached. Using either autodep or dependson will remove this warning. See the knitr cache options for more details.
#ADDENDUM: HRD signature analysis ```{r: HRD signature analysis, message=FALSE, warning=FALSE, echo=FALSE, cache=TRUE, eval=FALSE} sampleTable3 <- sampleTable %>% dplyr::filter(condition1==“Nigerian”) sampleTable3 <- sampleTable3 %>% dplyr::select(“sampleName”, “fileName”, “condition2”, “batch”) colnames(sampleTable3)[3] <- “sampleCondition”
ddsHTSeqMFNigerianonly <- DESeqDataSetFromHTSeqCount(sampleTable=sampleTable3, directory=FOLDER, design=~sampleCondition)
ddsHTSeqMFNigerianonly <- ddsHTSeqMFNigerianonly[rowSums(counts(ddsHTSeqMFNigerianonly)) > 0, ] #Pre-filtering the dataset by removing the rows without information about gene expression
dds <- estimateSizeFactors(ddsHTSeqMFNigerianonly)
HRDgenes=c(“ENSG00000149311.13”, “ENSG00000105393.11”, “ENSG00000138376.6”, “ENSG00000197299.6”, “ENSG00000012048.15”, “ENSG00000139618.10”, “ENSG00000158019.16”, “ENSG00000136492.4”, “ENSG00000183765.16”, “ENSG00000154920.10”, “ENSG00000163322.9”, “ENSG00000187741.10”, “ENSG00000187790.6”, “ENSG00000020922.8”, “ENSG00000104320.9”, “ENSG00000083093.5”, “ENSG00000062822.8”, “ENSG00000113522.9”, “ENSG00000182185.14”, “ENSG00000108384.10”, “ENSG00000185379.16”, “ENSG00000002016.12”, “ENSG00000197275.8”, “ENSG00000085999.7”, “ENSG00000101773.12”, “ENSG00000132383.7”, “ENSG00000139351.10”, “ENSG00000177302.10”, “ENSG00000100038.15”, “ENSG00000163781.8”, “ENSG00000087206.12”, “ENSG00000196584.2”, “ENSG00000126215.9”)
ATMcounts <- plotCounts(dds, gene=“ENSG00000149311.13”, intgroup=“sampleCondition”, returnData=TRUE) ATMcounts\(sampleName <- row.names(ATMcounts) ATMcounts\)percent10 <- ntile(ATMcounts$count, 10) temp <- ATMcounts %>% dplyr::filter(percent10==1) %>% dplyr::select(sampleName, sampleCondition, count) #write.table(temp, file=“ATMlowest10thpercentile.csv”, sep=“,”, row.names=FALSE, col.names=TRUE)
BARD1counts <- plotCounts(dds, gene=“ENSG00000138376.6”, intgroup=“sampleCondition”, returnData=TRUE) BARD1counts\(sampleName <- row.names(BARD1counts) BARD1counts\)percent10 <- ntile(BARD1counts$count, 10) temp <- BARD1counts %>% dplyr::filter(percent10==1) %>% dplyr::select(sampleName, sampleCondition, count) #write.table(temp, file=“BARD1lowest10thpercentile.csv”, sep=“,”, row.names=FALSE, col.names=TRUE)
BLMcounts <- plotCounts(dds, gene=“ENSG00000197299.6”, intgroup=“sampleCondition”, returnData=TRUE) BLMcounts\(sampleName <- row.names(BLMcounts) BLMcounts\)percent10 <- ntile(BLMcounts$count, 10) temp <- BLMcounts %>% dplyr::filter(percent10==1) %>% dplyr::select(sampleName, sampleCondition, count) #write.table(temp, file=“BLMlowest10thpercentile.csv”, sep=“,”, row.names=FALSE, col.names=TRUE)
BRCA1counts <- plotCounts(dds, gene=“ENSG00000012048.15”, intgroup=“sampleCondition”, returnData=TRUE) BRCA1counts\(sampleName <- row.names(BRCA1counts) BRCA1counts\)percent10 <- ntile(BRCA1counts$count, 10) temp <- BRCA1counts %>% dplyr::filter(percent10==1) %>% dplyr::select(sampleName, sampleCondition, count) #write.table(temp, file=“BRCA1lowest10thpercentile.csv”, sep=“,”, row.names=FALSE, col.names=TRUE)
BRCA2counts <- plotCounts(dds, gene=“ENSG00000139618.10”, intgroup=“sampleCondition”, returnData=TRUE) BRCA2counts\(sampleName <- row.names(BRCA2counts) BRCA2counts\)percent10 <- ntile(BRCA2counts$count, 10) temp <- BRCA2counts %>% dplyr::filter(percent10==1) %>% dplyr::select(sampleName, sampleCondition, count) #write.table(temp, file=“BRCA2lowest10thpercentile.csv”, sep=“,”, row.names=FALSE, col.names=TRUE)
BRIP1counts <- plotCounts(dds, gene=“ENSG00000136492.4”, intgroup=“sampleCondition”, returnData=TRUE) BRIP1counts\(sampleName <- row.names(BRIP1counts) BRIP1counts\)percent10 <- ntile(BRIP1counts$count, 10) temp <- BRIP1counts %>% dplyr::filter(percent10==1) %>% dplyr::select(sampleName, sampleCondition, count) #write.table(temp, file=“BRIP1lowest10thpercentile.csv”, sep=“,”, row.names=FALSE, col.names=TRUE)
CHEK2counts <- plotCounts(dds, gene=“ENSG00000183765.16”, intgroup=“sampleCondition”, returnData=TRUE) CHEK2counts\(sampleName <- row.names(CHEK2counts) CHEK2counts\)percent10 <- ntile(CHEK2counts$count, 10) temp <- CHEK2counts %>% dplyr::filter(percent10==1) %>% dplyr::select(sampleName, sampleCondition, count) #write.table(temp, file=“CHEK2lowest10thpercentile.csv”, sep=“,”, row.names=FALSE, col.names=TRUE)
FAM175Acounts <- plotCounts(dds, gene=“ENSG00000163322.9”, intgroup=“sampleCondition”, returnData=TRUE) FAM175Acounts\(sampleName <- row.names(FAM175Acounts) FAM175Acounts\)percent10 <- ntile(FAM175Acounts$count, 10) temp <- FAM175Acounts %>% dplyr::filter(percent10==1) %>% dplyr::select(sampleName, sampleCondition, count) #write.table(temp, file=“FAM175Alowest10thpercentile.csv”, sep=“,”, row.names=FALSE, col.names=TRUE)
FANCAcounts <- plotCounts(dds, gene=“ENSG00000187741.10”, intgroup=“sampleCondition”, returnData=TRUE) FANCAcounts\(sampleName <- row.names(FANCAcounts) FANCAcounts\)percent10 <- ntile(FANCAcounts$count, 10) temp <- FANCAcounts %>% dplyr::filter(percent10==1) %>% dplyr::select(sampleName, sampleCondition, count) #write.table(temp, file=“FANCAlowest10thpercentile.csv”, sep=“,”, row.names=FALSE, col.names=TRUE)
FANCMcounts <- plotCounts(dds, gene=“ENSG00000187790.6”, intgroup=“sampleCondition”, returnData=TRUE) FANCMcounts\(sampleName <- row.names(FANCMcounts) FANCMcounts\)percent10 <- ntile(FANCMcounts$count, 10) temp <- FANCMcounts %>% dplyr::filter(percent10==1) %>% dplyr::select(sampleName, sampleCondition, count) #write.table(temp, file=“FANCMlowest10thpercentile.csv”, sep=“,”, row.names=FALSE, col.names=TRUE)
MRE11Acounts <- plotCounts(dds, gene=“ENSG00000020922.8”, intgroup=“sampleCondition”, returnData=TRUE) MRE11Acounts\(sampleName <- row.names(MRE11Acounts) MRE11Acounts\)percent10 <- ntile(MRE11Acounts$count, 10) temp <- MRE11Acounts %>% dplyr::filter(percent10==1) %>% dplyr::select(sampleName, sampleCondition, count) #write.table(temp, file=“MRE11Alowest10thpercentile.csv”, sep=“,”, row.names=FALSE, col.names=TRUE)
NBNcounts <- plotCounts(dds, gene=“ENSG00000104320.9”, intgroup=“sampleCondition”, returnData=TRUE) NBNcounts\(sampleName <- row.names(NBNcounts) NBNcounts\)percent10 <- ntile(NBNcounts$count, 10) temp <- NBNcounts %>% dplyr::filter(percent10==1) %>% dplyr::select(sampleName, sampleCondition, count) #write.table(temp, file=“NBNlowest10thpercentile.csv”, sep=“,”, row.names=FALSE, col.names=TRUE)
PALB2counts <- plotCounts(dds, gene=“ENSG00000083093.5”, intgroup=“sampleCondition”, returnData=TRUE) PALB2counts\(sampleName <- row.names(PALB2counts) PALB2counts\)percent10 <- ntile(PALB2counts$count, 10) temp <- PALB2counts %>% dplyr::filter(percent10==1) %>% dplyr::select(sampleName, sampleCondition, count) #write.table(temp, file=“PALB2lowest10thpercentile.csv”, sep=“,”, row.names=FALSE, col.names=TRUE)
POLD1counts <- plotCounts(dds, gene=“ENSG00000062822.8”, intgroup=“sampleCondition”, returnData=TRUE) POLD1counts\(sampleName <- row.names(POLD1counts) POLD1counts\)percent10 <- ntile(POLD1counts$count, 10) temp <- POLD1counts %>% dplyr::filter(percent10==1) %>% dplyr::select(sampleName, sampleCondition, count) #write.table(temp, file=“POLD1lowest10thpercentile.csv”, sep=“,”, row.names=FALSE, col.names=TRUE)
RAD51Ccounts <- plotCounts(dds, gene=“ENSG00000108384.10”, intgroup=“sampleCondition”, returnData=TRUE) RAD51Ccounts\(sampleName <- row.names(RAD51Ccounts) RAD51Ccounts\)percent10 <- ntile(RAD51Ccounts$count, 10) temp <- RAD51Ccounts %>% dplyr::filter(percent10==1) %>% dplyr::select(sampleName, sampleCondition, count) #write.table(temp, file=“RAD51Clowest10thpercentile.csv”, sep=“,”, row.names=FALSE, col.names=TRUE)
XRCCcounts <- plotCounts(dds, gene=“ENSG00000196584.2”, intgroup=“sampleCondition”, returnData=TRUE) XRCCcounts\(sampleName <- row.names(XRCCcounts) XRCCcounts\)percent10 <- ntile(XRCCcounts$count, 10) temp <- XRCCcounts %>% dplyr::filter(percent10==1) %>% dplyr::select(sampleName, sampleCondition, count) #write.table(temp, file=“XRCClowest10thpercentile.csv”, sep=“,”, row.names=FALSE, col.names=TRUE)
testgenename=“ENSG00000149311.13” plotCounts(dds, gene = testgenename, intgroup=“sampleCondition”, main=“Distribution of ATM expression in Nigerian cohort”)
testgenename=“ENSG00000105393.11” plotCounts(dds, gene = testgenename, intgroup=“sampleCondition”, main=“Distribution of BABAM1 expression in Nigerian cohort”)
testgenename=“ENSG00000138376.6” plotCounts(dds, gene = testgenename, intgroup=“sampleCondition”, main=“Distribution of BARD1 expression in Nigerian cohort”)
testgenename=“ENSG00000197299.6” plotCounts(dds, gene = testgenename, intgroup=“sampleCondition”, main=“Distribution of BLM expression in Nigerian cohort”)
testgenename=“ENSG00000012048.15” plotCounts(dds, gene = testgenename, intgroup=“sampleCondition”, main=“Distribution of BRCA1 expression in Nigerian cohort”)
testgenename=“ENSG00000139618.10” plotCounts(dds, gene = testgenename, intgroup=“sampleCondition”, main=“Distribution of BRCA2 expression in Nigerian cohort”)
testgenename=“ENSG00000158019.16” plotCounts(dds, gene = testgenename, intgroup=“sampleCondition”, main=“Distribution of BRE expression in Nigerian cohort”)
testgenename=“ENSG00000136492.4” plotCounts(dds, gene = testgenename, intgroup=“sampleCondition”, main=“Distribution of BRIP1 expression in Nigerian cohort”)
testgenename=“ENSG00000183765.16” plotCounts(dds, gene = testgenename, intgroup=“sampleCondition”, main=“Distribution of CHEK2 expression in Nigerian cohort”)
testgenename=“ENSG00000154920.10” plotCounts(dds, gene = testgenename, intgroup=“sampleCondition”, main=“Distribution of EME1 expression in Nigerian cohort”)
testgenename=“ENSG00000163322.9” plotCounts(dds, gene = testgenename, intgroup=“sampleCondition”, main=“Distribution of FAM175A expression in Nigerian cohort”)
testgenename=“ENSG00000187741.10” plotCounts(dds, gene = testgenename, intgroup=“sampleCondition”, main=“Distribution of FANCA expression in Nigerian cohort”)
testgenename=“ENSG00000187790.6” plotCounts(dds, gene = testgenename, intgroup=“sampleCondition”, main=“Distribution of FANCM expression in Nigerian cohort”)
testgenename=“ENSG00000020922.8” plotCounts(dds, gene = testgenename, intgroup=“sampleCondition”, main=“Distribution of MRE11A expression in Nigerian cohort”)
testgenename=“ENSG00000104320.9” plotCounts(dds, gene = testgenename, intgroup=“sampleCondition”, main=“Distribution of NBN expression in Nigerian cohort”)
testgenename=“ENSG00000083093.5” plotCounts(dds, gene = testgenename, intgroup=“sampleCondition”, main=“Distribution of PALB2 expression in Nigerian cohort”)
testgenename=“ENSG00000062822.8” plotCounts(dds, gene = testgenename, intgroup=“sampleCondition”, main=“Distribution of POLD1 expression in Nigerian cohort”)
testgenename=“ENSG00000113522.9” plotCounts(dds, gene = testgenename, intgroup=“sampleCondition”, main=“Distribution of RAD50 expression in Nigerian cohort”)
testgenename=“ENSG00000182185.14” plotCounts(dds, gene = testgenename, intgroup=“sampleCondition”, main=“Distribution of RAD51B expression in Nigerian cohort”)
testgenename=“ENSG00000108384.10” plotCounts(dds, gene = testgenename, intgroup=“sampleCondition”, main=“Distribution of RAD51C expression in Nigerian cohort”)
testgenename=“ENSG00000185379.16” plotCounts(dds, gene = testgenename, intgroup=“sampleCondition”, main=“Distribution of RAD51D expression in Nigerian cohort”)
testgenename=“ENSG00000002016.12” plotCounts(dds, gene = testgenename, intgroup=“sampleCondition”, main=“Distribution of RAD52 expression in Nigerian cohort”)
testgenename=“ENSG00000197275.8” plotCounts(dds, gene = testgenename, intgroup=“sampleCondition”, main=“Distribution of RAD54B expression in Nigerian cohort”)
testgenename=“ENSG00000085999.7” plotCounts(dds, gene = testgenename, intgroup=“sampleCondition”, main=“Distribution of RAD54L expression in Nigerian cohort”)
testgenename=“ENSG00000101773.12” plotCounts(dds, gene = testgenename, intgroup=“sampleCondition”, main=“Distribution of RBBP8 expression in Nigerian cohort”)
testgenename=“ENSG00000132383.7” plotCounts(dds, gene = testgenename, intgroup=“sampleCondition”, main=“Distribution of RPA1 expression in Nigerian cohort”)
testgenename=“ENSG00000139351.10” plotCounts(dds, gene = testgenename, intgroup=“sampleCondition”, main=“Distribution of SYCP3 expression in Nigerian cohort”)
testgenename=“ENSG00000177302.10” plotCounts(dds, gene = testgenename, intgroup=“sampleCondition”, main=“Distribution of TOP3A expression in Nigerian cohort”)
testgenename=“ENSG00000100038.15” plotCounts(dds, gene = testgenename, intgroup=“sampleCondition”, main=“Distribution of TOP3B expression in Nigerian cohort”)
testgenename=“ENSG00000163781.8” plotCounts(dds, gene = testgenename, intgroup=“sampleCondition”, main=“Distribution of TOPBP1 expression in Nigerian cohort”)
testgenename=“ENSG00000087206.12” plotCounts(dds, gene = testgenename, intgroup=“sampleCondition”, main=“Distribution of UIMC1 expression in Nigerian cohort”)
testgenename=“ENSG00000196584.2” plotCounts(dds, gene = testgenename, intgroup=“sampleCondition”, main=“Distribution of XRCC2 expression in Nigerian cohort”)
testgenename=“ENSG00000126215.9” plotCounts(dds, gene = testgenename, intgroup=“sampleCondition”, main=“Distribution of XRCC3 expression in Nigerian cohort”) ```
sessionInfo()R version 3.6.3 (2020-02-29)
Platform: x86_64-apple-darwin15.6.0 (64-bit)
Running under: macOS Sierra 10.12.6
Matrix products: default
BLAS: /Library/Frameworks/R.framework/Versions/3.6/Resources/lib/libRblas.0.dylib
LAPACK: /Library/Frameworks/R.framework/Versions/3.6/Resources/lib/libRlapack.dylib
locale:
[1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
attached base packages:
[1] parallel stats4 grid stats graphics grDevices utils
[8] datasets methods base
other attached packages:
[1] Glimma_1.12.0 RColorBrewer_1.1-2
[3] preprocessCore_1.46.0 ashr_2.2-32
[5] ggfortify_0.4.7 calibrate_1.7.2
[7] MASS_7.3-51.5 sva_3.32.1
[9] mgcv_1.8-31 nlme_3.1-144
[11] EnsDb.Hsapiens.v75_2.99.0 ensembldb_2.8.0
[13] AnnotationFilter_1.8.0 GenomicFeatures_1.36.4
[15] hexbin_1.27.3 stringi_1.4.3
[17] dplyr_0.8.3 affy_1.62.0
[19] checkmate_1.9.3 pathview_1.24.0
[21] org.Hs.eg.db_3.8.2 AnnotationDbi_1.46.0
[23] clusterProfiler_3.12.0 pheatmap_1.0.12
[25] genefilter_1.66.0 vsn_3.52.0
[27] RUVSeq_1.18.0 EDASeq_2.18.0
[29] ShortRead_1.42.0 GenomicAlignments_1.20.0
[31] Rsamtools_2.0.0 Biostrings_2.52.0
[33] XVector_0.24.0 DESeq2_1.24.0
[35] SummarizedExperiment_1.14.0 DelayedArray_0.10.0
[37] BiocParallel_1.18.0 matrixStats_0.54.0
[39] Biobase_2.44.0 GenomicRanges_1.36.0
[41] GenomeInfoDb_1.20.0 IRanges_2.18.1
[43] S4Vectors_0.22.0 BiocGenerics_0.30.0
[45] edgeR_3.26.4 limma_3.40.2
[47] ggbiplot_0.55 scales_1.0.0
[49] plyr_1.8.5 ggplot2_3.2.1
[51] gplots_3.0.3
loaded via a namespace (and not attached):
[1] R.utils_2.8.0 tidyselect_0.2.5 RSQLite_2.1.1
[4] htmlwidgets_1.3 DESeq_1.36.0 munsell_0.5.0
[7] codetools_0.2-16 withr_2.1.2 colorspace_1.4-1
[10] GOSemSim_2.10.0 knitr_1.28 rstudioapi_0.11
[13] pscl_1.5.2 DOSE_3.10.1 labeling_0.3
[16] git2r_0.26.1 KEGGgraph_1.44.0 urltools_1.7.3
[19] GenomeInfoDbData_1.2.1 mixsqp_0.1-97 hwriter_1.3.2
[22] polyclip_1.10-0 bit64_0.9-7 farver_1.1.0
[25] rprojroot_1.3-2 vctrs_0.2.0 xfun_0.7
[28] doParallel_1.0.14 R6_2.4.0 locfit_1.5-9.1
[31] bitops_1.0-6 fgsea_1.10.0 gridGraphics_0.4-1
[34] assertthat_0.2.1 ggraph_1.0.2 nnet_7.3-12
[37] enrichplot_1.4.0 gtable_0.3.0 workflowr_1.4.0
[40] rlang_0.4.5 zeallot_0.1.0 splines_3.6.3
[43] rtracklayer_1.44.0 lazyeval_0.2.2 acepack_1.4.1
[46] europepmc_0.3 BiocManager_1.30.10 yaml_2.2.0
[49] reshape2_1.4.3 backports_1.1.4 qvalue_2.16.0
[52] Hmisc_4.2-0 tools_3.6.3 ggplotify_0.0.3
[55] affyio_1.54.0 ggridges_0.5.1 Rcpp_1.0.1
[58] base64enc_0.1-3 progress_1.2.2 zlibbioc_1.30.0
[61] purrr_0.3.3 RCurl_1.95-4.12 prettyunits_1.0.2
[64] rpart_4.1-15 viridis_0.5.1 cowplot_0.9.4
[67] ggrepel_0.8.1 cluster_2.1.0 fs_1.3.1
[70] magrittr_1.5 data.table_1.12.8 DO.db_2.9
[73] triebeard_0.3.0 truncnorm_1.0-8 SQUAREM_2017.10-1
[76] whisker_0.3-2 ProtGenerics_1.16.0 aroma.light_3.14.0
[79] hms_0.5.2 evaluate_0.14 xtable_1.8-4
[82] XML_3.98-1.20 gridExtra_2.3 compiler_3.6.3
[85] biomaRt_2.40.0 tibble_2.1.3 KernSmooth_2.23-16
[88] crayon_1.3.4 R.oo_1.22.0 htmltools_0.3.6
[91] Formula_1.2-3 tidyr_1.0.0 geneplotter_1.62.0
[94] DBI_1.0.0 tweenr_1.0.1 Matrix_1.2-18
[97] R.methodsS3_1.7.1 gdata_2.18.0 igraph_1.2.4.1
[100] pkgconfig_2.0.2 rvcheck_0.1.3 foreign_0.8-75
[103] foreach_1.4.4 xml2_1.3.2 annotate_1.62.0
[106] stringr_1.4.0 digest_0.6.25 graph_1.62.0
[109] rmarkdown_2.1 fastmatch_1.1-0 htmlTable_1.13.1
[112] curl_4.3 gtools_3.8.1 lifecycle_0.1.0
[115] jsonlite_1.6.1 viridisLite_0.3.0 pillar_1.4.2
[118] lattice_0.20-38 KEGGREST_1.24.0 httr_1.4.1
[121] survival_3.1-8 GO.db_3.8.2 glue_1.4.0
[124] UpSetR_1.4.0 iterators_1.0.10 png_0.1-7
[127] bit_1.1-14 Rgraphviz_2.28.0 ggforce_0.2.2
[130] blob_1.1.1 latticeExtra_0.6-28 caTools_1.17.1.2
[133] memoise_1.1.0